% PASADENA experiment simulation for the parahydrogenation of styrene
% into ethylbenzene.
%
% i.kuprov@soton.ac.uk

function pasadena_ethylbenzene()

% Spin system
sys.isotopes={'1H','1H','1H','1H','1H','1H','1H','1H','1H','1H'};

% Magnetic field
sys.magnet=7.05;

% Chemical shifts
inter.zeeman.scalar={1.201 1.201 1.201 2.625 2.625...
                     7.207 7.207 7.265 7.265 7.155};

% Scalar couplings
inter.coupling.scalar=cell(10);
inter.coupling.scalar{1,4}=7.63;   inter.coupling.scalar{2,4}=7.63;
inter.coupling.scalar{3,4}=7.63;   inter.coupling.scalar{1,5}=7.63;
inter.coupling.scalar{2,5}=7.63;   inter.coupling.scalar{3,5}=7.63;
inter.coupling.scalar{4,6}=-0.62;  inter.coupling.scalar{4,7}=-0.62;
inter.coupling.scalar{5,6}=-0.62;  inter.coupling.scalar{5,7}=-0.62;
inter.coupling.scalar{4,8}=0.29;   inter.coupling.scalar{4,9}=0.29;
inter.coupling.scalar{5,8}=0.29;   inter.coupling.scalar{5,9}=0.29;
inter.coupling.scalar{4,10}=-0.45; inter.coupling.scalar{5,10}=-0.45;
inter.coupling.scalar{6,8}=7.63;   inter.coupling.scalar{6,10}=1.26;
inter.coupling.scalar{7,9}=7.63;   inter.coupling.scalar{7,10}=1.26;
inter.coupling.scalar{8,10}=7.44;  inter.coupling.scalar{9,10}=7.44;

% Basis set
bas.mode='complete';
bas.sym_group={'S3','S2','S2','S2'};
bas.sym_spins={[1 2 3],[4 5],[6 7],[8 9]};

% Spinach housekeeping
spin_system=create(sys,inter);
spin_system=basis(spin_system,bas);

% Sequence parameters
parameters.spins={'1H'};
parameters.rho_init=state(spin_system,{'Lz','Lz'},{1,4});
parameters.coil=state(spin_system,'L+','1H');
parameters.pulse_op=(operator(spin_system,'L+','1H')-operator(spin_system,'L-','1H'))/2i;
parameters.pulse_angle=pi/4;
parameters.offset=500;
parameters.sweep=1000;
parameters.npoints=1024;
parameters.zerofill=8192;
parameters.axis_units='ppm';
parameters.invert_axis=1;

% Simulation
fid=liquid(spin_system,@pulse_acquire,parameters,'nmr');

% Apodization
fid=apodization(fid,'gaussian-1d',10);

% Fourier transform
spectrum=fftshift(fft(fid,parameters.zerofill));

% Plotting
plot_1d(spin_system,real(spectrum),parameters);

end

